{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "940844d8",
   "metadata": {},
   "source": [
    "在minist手写数字上训练hvae图像生成模型"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5473f49f",
   "metadata": {},
   "outputs": [],
   "source": [
    "import torch\n",
    "import torchvision as tv\n",
    "import torch.nn as nn\n",
    "import torch.nn.functional as F\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "2e2b6cf6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# 超参数定义\n",
    "from dataclasses import dataclass\n",
    "\n",
    "\n",
    "@dataclass\n",
    "class Config:\n",
    "    # 文件&数据\n",
    "    data_path: str = \"./dataset/mnist\"\n",
    "    save_dir: str = \"./checkpoints\"\n",
    "    save_every: int = 10\n",
    "    transform = tv.transforms.Compose([\n",
    "        tv.transforms.ToTensor(),\n",
    "        # tv.transforms.Normalize((0.5,), (0.5,))\n",
    "    ])\n",
    "    # 训练参数\n",
    "    device: str = \"cuda\" if torch.cuda.is_available() else \"cpu\"\n",
    "    batch_size: int = 1024\n",
    "    num_workers: int = 4\n",
    "    lr: float = 1e-4\n",
    "    num_epochs: int = 300\n",
    "    # 模型参数\n",
    "    img_size = 28\n",
    "    latent_dim1 = 32\n",
    "    latent_dim2 = 32\n",
    "\n",
    "g_opt = Config()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "36bf4028",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "共有60000张训练图片，10000张测试图片\n"
     ]
    }
   ],
   "source": [
    "# 数据集\n",
    "from torch.utils.data import DataLoader\n",
    "\n",
    "\n",
    "train_dataset = tv.datasets.MNIST(\n",
    "    root=g_opt.data_path,\n",
    "    train=True,\n",
    "    transform=g_opt.transform,\n",
    "    download=True\n",
    ")\n",
    "test_dataset = tv.datasets.MNIST(\n",
    "    root=g_opt.data_path,\n",
    "    train=False,\n",
    "    transform=g_opt.transform,\n",
    "    download=True\n",
    ")\n",
    "print(f\"共有{len(train_dataset)}张训练图片，{len(test_dataset)}张测试图片\")\n",
    "# 数据加载器\n",
    "train_loader = DataLoader(\n",
    "    train_dataset,\n",
    "    batch_size=g_opt.batch_size,\n",
    "    shuffle=True,\n",
    "    num_workers=g_opt.num_workers\n",
    ")\n",
    "val_loader = DataLoader(\n",
    "    test_dataset,\n",
    "    batch_size=g_opt.batch_size,\n",
    "    shuffle=False,\n",
    "    num_workers=g_opt.num_workers\n",
    ")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a6c52e13",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Utility\n",
    "def save_images(tensor, fn, nrow=8, padding=2):\n",
    "    \"\"\"\n",
    "    保存图像\n",
    "    \"\"\"\n",
    "    grid = tv.utils.make_grid(tensor, nrow=nrow, padding=padding)\n",
    "    tv.utils.save_image(grid, fn)\n",
    "\n",
    "def kl_divergence_gaussians(q_mu, q_logvar, p_mu, p_logvar):\n",
    "    \"\"\"\n",
    "    计算两个多元高斯分布之间的KL散度 D_KL(q || p)\n",
    "    q: N(q_mu, exp(q_logvar))\n",
    "    p: N(p_mu, exp(p_logvar))\n",
    "    \"\"\"\n",
    "    q_std = torch.exp(0.5 * q_logvar)\n",
    "    p_std = torch.exp(0.5 * p_logvar)\n",
    "    \n",
    "    # D_KL(q || p) = 0.5 * [ log(det(Sigma_p)/det(Sigma_q)) - k + tr(Sigma_p^-1 * Sigma_q) + (mu_p - mu_q)^T Sigma_p^-1 (mu_p - mu_q) ]\n",
    "    # 对于对角协方差矩阵，公式简化为：\n",
    "    kl_div = torch.sum(p_logvar - q_logvar - 1 + (q_std.pow(2) + (q_mu - p_mu).pow(2)) / p_std.pow(2), dim=1)\n",
    "    \n",
    "    return 0.5 * kl_div\n",
    "\n",
    "def hvae_loss(x, x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1):\n",
    "    # 1. 重构损失\n",
    "    recon_loss = F.binary_cross_entropy(x_recon.view(-1, 784), x.view(-1, 784), reduction='sum')\n",
    "    # 2. KL 散度 for z2\n",
    "    # D_KL(q(z2|z1) || p(z2)), 其中 p(z2) ~ N(0, I)\n",
    "    mu_p2 = torch.zeros_like(mu_q2)\n",
    "    logvar_p2 = torch.zeros_like(logvar_q2)\n",
    "    kl2 = kl_divergence_gaussians(mu_q2, logvar_q2, mu_p2, logvar_p2).sum()\n",
    "    # 3. KL 散度 for z1\n",
    "    # D_KL(q(z1|x) || p(z1|z2))\n",
    "    kl1 = kl_divergence_gaussians(mu_q1, logvar_q1, mu_p1, logvar_p1).sum()\n",
    "    batch_size = x.size(0)\n",
    "    kl1 /= batch_size\n",
    "    kl2 /= batch_size\n",
    "    recon_loss /= batch_size\n",
    "    return recon_loss + kl1 + kl2, recon_loss, kl1 + kl2\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "81cbece5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# HVAE Model\n",
    "from typing import Tuple\n",
    "\n",
    "\n",
    "class HVAE(nn.Module):\n",
    "    def __init__(self, img_size=28, latent_dim1=20, latent_dim2=10):\n",
    "        super(HVAE, self).__init__()\n",
    "        self.img_size = img_size\n",
    "        self.latent_dim1 = latent_dim1\n",
    "        self.latent_dim2 = latent_dim2\n",
    "        \n",
    "        # --- 1. 推断模型 (Inference Model: q(z1|x), q(z2|z1)) ---\n",
    "        \n",
    "        # 1a. 从 x 提取特征\n",
    "        self.encoder_base = nn.Sequential(\n",
    "            nn.Conv2d(1, 16, kernel_size=4, stride=2, padding=1), # (B, 1, 28, 28) -> (B, 16, 14, 14)\n",
    "            nn.ReLU(),\n",
    "            nn.Conv2d(16, 32, kernel_size=4, stride=2, padding=1), # (B, 16, 14, 14) -> (B, 32, 7, 7)\n",
    "            nn.ReLU(),\n",
    "            nn.Conv2d(32, 64, kernel_size=3, stride=2, padding=1), # (B, 32, 7, 7) -> (B, 64, 4, 4)\n",
    "            nn.ReLU(),\n",
    "            nn.Flatten() # (B, 64*4*4=1024)\n",
    "        )\n",
    "        \n",
    "        # 1b. 从特征推断 z1: q(z1|x)\n",
    "        self.fc_q1_mu = nn.Linear(64 * 4 * 4, latent_dim1)\n",
    "        self.fc_q1_logvar = nn.Linear(64 * 4 * 4, latent_dim1)\n",
    "        \n",
    "        # 1c. 从 z1 推断 z2: q(z2|z1)\n",
    "        self.encoder_z2 = nn.Sequential(\n",
    "            nn.Linear(latent_dim1, 128),\n",
    "            nn.ReLU()\n",
    "        )\n",
    "        self.fc_q2_mu = nn.Linear(128, latent_dim2)\n",
    "        self.fc_q2_logvar = nn.Linear(128, latent_dim2)\n",
    "        # --- 2. 生成模型 (Generative Model: p(z1|z2), p(x|z1)) ---\n",
    "        # 2a. 从 z2 生成 z1 的先验分布: p(z1|z2)\n",
    "        self.decoder_p1 = nn.Sequential(\n",
    "            nn.Linear(latent_dim2, 128),\n",
    "            nn.ReLU()\n",
    "        )\n",
    "        self.fc_p1_mu = nn.Linear(128, latent_dim1)\n",
    "        self.fc_p1_logvar = nn.Linear(128, latent_dim1)\n",
    "        \n",
    "        # 2b. 从 z1 重构 x: p(x|z1)\n",
    "        self.decoder_input = nn.Linear(latent_dim1, 64 * 4 * 4)\n",
    "        self.decoder_x = nn.Sequential(\n",
    "            nn.Unflatten(1, (64, 4, 4)),\n",
    "            nn.ConvTranspose2d(64, 32, kernel_size=3, stride=2, padding=1), # (B, 64, 4, 4) -> (B, 32, 7, 7)\n",
    "            nn.ReLU(),\n",
    "            nn.ConvTranspose2d(32, 16, kernel_size=4, stride=2, padding=1), # (B, 32, 7, 7) -> (B, 16, 14, 14)\n",
    "            nn.ReLU(),\n",
    "            nn.ConvTranspose2d(16, 1, kernel_size=4, stride=2, padding=1), # (B, 16, 14, 14) -> (B, 1, 28, 28)\n",
    "            nn.Sigmoid(),\n",
    "        )\n",
    "    def reparameterize(self, mu, logvar):\n",
    "        std = torch.exp(0.5 * logvar)\n",
    "        eps = torch.randn_like(std)\n",
    "        return mu + eps * std\n",
    "    def encode(self, x: torch.Tensor) -> Tuple[torch.Tensor, ...]:\n",
    "        # 推断 q(z1|x) 和 q(z2|z1)\n",
    "        h = self.encoder_base(x)\n",
    "        \n",
    "        # q(z1|x)\n",
    "        mu_q1 = self.fc_q1_mu(h)\n",
    "        logvar_q1 = self.fc_q1_logvar(h)\n",
    "        z1 = self.reparameterize(mu_q1, logvar_q1)\n",
    "        \n",
    "        # q(z2|z1)\n",
    "        h2 = self.encoder_z2(z1)\n",
    "        mu_q2 = self.fc_q2_mu(h2)\n",
    "        logvar_q2 = self.fc_q2_logvar(h2)\n",
    "        z2 = self.reparameterize(mu_q2, logvar_q2)\n",
    "        \n",
    "        return z1, mu_q1, logvar_q1, z2, mu_q2, logvar_q2\n",
    "    def decode(self, z1: torch.Tensor, z2: torch.Tensor) -> Tuple[torch.Tensor, ...]:\n",
    "        # 生成 p(z1|z2) 和 p(x|z1)\n",
    "        \n",
    "        # p(z1|z2)\n",
    "        h_p1 = self.decoder_p1(z2)\n",
    "        mu_p1 = self.fc_p1_mu(h_p1)\n",
    "        logvar_p1 = self.fc_p1_logvar(h_p1)\n",
    "        \n",
    "        # p(x|z1)\n",
    "        h_x = self.decoder_input(z1)\n",
    "        x_recon = self.decoder_x(h_x)\n",
    "        \n",
    "        return x_recon, mu_p1, logvar_p1\n",
    "    def forward(self, x: torch.Tensor):\n",
    "        # 推断路径\n",
    "        z1, mu_q1, logvar_q1, z2, mu_q2, logvar_q2 = self.encode(x)\n",
    "        \n",
    "        # 生成路径\n",
    "        x_recon, mu_p1, logvar_p1 = self.decode(z1, z2)\n",
    "        \n",
    "        return x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1\n",
    "        \n",
    "    def sample(self, num_samples: int, device: torch.device):\n",
    "        \"\"\"从先验中采样生成新图像\"\"\"\n",
    "        # 从顶层先验 p(z2) ~ N(0, I) 采样\n",
    "        z2 = torch.randn(num_samples, self.latent_dim2).to(device)\n",
    "        \n",
    "        # 基于 z2 计算 p(z1|z2) 的参数\n",
    "        h_p1 = self.decoder_p1(z2)\n",
    "        mu_p1 = self.fc_p1_mu(h_p1)\n",
    "        logvar_p1 = self.fc_p1_logvar(h_p1)\n",
    "        \n",
    "        # 从 p(z1|z2) 采样 z1\n",
    "        z1 = self.reparameterize(mu_p1, logvar_p1)\n",
    "        \n",
    "        # 从 z1 生成图像\n",
    "        h_x = self.decoder_input(z1)\n",
    "        generated_x = self.decoder_x(h_x)\n",
    "        \n",
    "        return generated_x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "21652e20",
   "metadata": {},
   "outputs": [],
   "source": [
    "from tqdm import tqdm\n",
    "\n",
    "def train_epoch(model, train_loader, optimizer, device):\n",
    "    model.train()\n",
    "    train_loss, recon_loss, kld_loss = 0, 0, 0\n",
    "    for x, y in train_loader:\n",
    "        x = x.to(device)\n",
    "        x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1 = model(x)\n",
    "        loss, recon_loss, kld_loss = hvae_loss(x, x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1)\n",
    "        optimizer.zero_grad()\n",
    "        loss.backward()\n",
    "        optimizer.step()\n",
    "        train_loss += loss.item()\n",
    "        recon_loss += recon_loss.item()\n",
    "        kld_loss += kld_loss.item()\n",
    "        \n",
    "    return train_loss/len(train_loader), recon_loss/len(train_loader), kld_loss/len(train_loader)\n",
    "\n",
    "def val_epoch(model, val_loader, device):\n",
    "    model.eval()\n",
    "    with torch.no_grad():\n",
    "        val_loss, recon_loss, kld_loss = 0, 0, 0\n",
    "        for x, y in val_loader:\n",
    "            x = x.to(device)\n",
    "            x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1 = model(x)\n",
    "            loss, recon_loss, kld_loss = hvae_loss(x, x_recon, mu_q1, logvar_q1, mu_q2, logvar_q2, mu_p1, logvar_p1)\n",
    "            val_loss += loss.item()\n",
    "            recon_loss += recon_loss.item()\n",
    "            kld_loss += kld_loss.item()\n",
    "            \n",
    "    return val_loss/len(val_loader), recon_loss/len(val_loader), kld_loss/len(val_loader)\n",
    "\n",
    "def train(model, train_loader, val_loader, optimizer, device, num_epochs):\n",
    "    model.to(device)\n",
    "    best_loss = float('inf')\n",
    "\n",
    "    epoch_pbar = tqdm(range(num_epochs), desc=\"Epochs\", position=0)\n",
    "\n",
    "    for epoch in epoch_pbar:\n",
    "        train_loss, _, _ = train_epoch(model, train_loader, optimizer, device)\n",
    "        val_loss, recon_loss, kld_loss = val_epoch(model, val_loader, device)\n",
    "        epoch_pbar.set_postfix(\n",
    "            Train=f\"{train_loss:.4f}\",\n",
    "            Test=f\"{val_loss:.4f}\",\n",
    "            Recon=f\"{recon_loss:.3f}\",\n",
    "            KLD=f\"{kld_loss:.3f}\"\n",
    "        )\n",
    "        # print(f\"Epoch {epoch+1}/{num_epochs} TrainLoss: {train_loss:.4f} TestLoss: {val_loss:.4f} ReconLoss: {recon_loss:.3f} KLD: {kld_loss:.3f}\")        \n",
    "        cur_val_loss = val_loss\n",
    "        # 保存最佳模型\n",
    "        if cur_val_loss < best_loss:\n",
    "            best_loss = cur_val_loss\n",
    "            torch.save({\n",
    "                'epoch': epoch,\n",
    "                'model_state_dict': model.state_dict(),\n",
    "                'optimizer_state_dict': optimizer.state_dict(),\n",
    "                'loss': best_loss\n",
    "            }, f\"{g_opt.save_dir}/best.pth\")\n",
    "            \n",
    "        torch.save({\n",
    "            'epoch': epoch,\n",
    "            'model_state_dict': model.state_dict(),\n",
    "            'optimizer_state_dict': optimizer.state_dict(),\n",
    "            'loss': cur_val_loss\n",
    "        }, f\"{g_opt.save_dir}/last.pth\")\n",
    "\n",
    "        # if epoch % g_opt.save_every == 0 or epoch == num_epochs - 1 or epoch == 0:\n",
    "        #     save_epoch(model, val_loader, epoch, device, times=1)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "f3db9372",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Epochs: 100%|██████████| 300/300 [03:24<00:00,  1.47it/s, KLD=4040.264, Recon=11891.041, Test=100400.2961, Train=102824.1690]\n"
     ]
    }
   ],
   "source": [
    "model = HVAE(img_size=g_opt.img_size, latent_dim1=g_opt.latent_dim1, latent_dim2=g_opt.latent_dim2)\n",
    "optimizer = torch.optim.Adam(model.parameters(), lr=g_opt.lr)\n",
    "train(model, train_loader, val_loader, optimizer, g_opt.device, g_opt.num_epochs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e270b081",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "测试集图像重建效果：\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7YAAAE1CAYAAADTQXzHAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAATlBJREFUeJzt3Xd8VMX6+PEnoYQWQknoPaFFkM5FOoiiAgIKKiISQJpS/F5AUS5VBUFRvIiAFRRFpYMogkoRwQJIEUSa9JqEEjok5/eHP+bOTJIlhN0kZ/N5v16+Xs/sc7I77Nlzdsczz5kAx3EcAQAAAADApQLTuwMAAAAAANwOBrYAAAAAAFdjYAsAAAAAcDUGtgAAAAAAV2NgCwAAAABwNQa2AAAAAABXY2ALAAAAAHA1BrYAAAAAAFdjYAsAAAAAcDUGtgAAZGCjRo2SgICAVP3tjBkzJCAgQPbv3+/dTmn2798vAQEBMmPGDJ+9BgAAN8PAFgAAH9i+fbs88cQTUrx4cQkKCpJixYpJ586dZfv27endtXSxatUqCQgIkLlz56Z3VwAAfoiBLQAAXjZ//nypWbOmfP/999KtWzd55513pEePHrJy5UqpWbOmLFiwIMXP9Z///EcuXbqUqn506dJFLl26JKVLl07V3wMA4BZZ07sDAAD4k71790qXLl2kXLlysmbNGgkLC1O5gQMHSqNGjaRLly6ydetWKVeuXLLPc+HCBcmdO7dkzZpVsmZN3dd1lixZJEuWLKn6WwAA3IQrtgAAeNFrr70mFy9elHfffdcY1IqIhIaGyvTp0+XChQsyYcIE9fiNOtodO3bI448/Lvnz55eGDRsaOd2lS5dkwIABEhoaKsHBwfLggw/KkSNHJCAgQEaNGqW2S6rGtkyZMtK6dWtZu3at1K1bV3LkyCHlypWTjz/+2HiN2NhYGTx4sFStWlXy5MkjefPmlfvvv1+2bNnipXfqf/+2Xbt2yRNPPCEhISESFhYmw4cPF8dx5NChQ9K2bVvJmzevFClSRCZOnGj8/dWrV2XEiBFSq1YtCQkJkdy5c0ujRo1k5cqViV4rJiZGunTpInnz5pV8+fJJ165dZcuWLUnWB+/cuVM6dOggBQoUkBw5ckjt2rVl8eLFxjbXrl2T0aNHS/ny5SVHjhxSsGBBadiwoaxYscJr7w8AIOUY2AIA4EVLliyRMmXKSKNGjZLMN27cWMqUKSNLly5NlOvYsaNcvHhRxo4dKz179kz2NaKiomTy5MnywAMPyPjx4yVnzpzSqlWrFPdxz5490qFDB7nnnntk4sSJkj9/fomKijLqf/ft2ycLFy6U1q1byxtvvCFDhgyRbdu2SZMmTeTo0aMpfq2UePTRRyUhIUFeffVV+de//iUvv/yyTJo0Se655x4pXry4jB8/XiIiImTw4MGyZs0a9Xfnzp2T999/X5o2bSrjx4+XUaNGyalTp6Rly5ayefNmtV1CQoK0adNGZs+eLV27dpVXXnlFjh07Jl27dk3Ul+3bt0u9evXkzz//lKFDh8rEiRMld+7c0q5dO2MK+ahRo2T06NHSrFkzefvtt2XYsGFSqlQp2bRpk1ffGwBACjkAAMArzpw544iI07ZtW4/bPfjgg46IOOfOnXMcx3FGjhzpiIjTqVOnRNveyN2wceNGR0ScZ5991tguKirKERFn5MiR6rGPPvrIERHn77//Vo+VLl3aERFnzZo16rGTJ086QUFBzqBBg9Rjly9fduLj443X+Pvvv52goCBnzJgxxmMi4nz00Uce/80rV650RMSZM2dOon9br1691GPXr193SpQo4QQEBDivvvqqevz06dNOzpw5na5duxrbXrlyxXid06dPO4ULF3a6d++uHps3b54jIs6kSZPUY/Hx8U7z5s0T9f3uu+92qlat6ly+fFk9lpCQ4NSvX98pX768eqxatWpOq1atPP6bAQBphyu2AAB4SVxcnIiIBAcHe9zuRv7cuXPG43369LnpayxbtkxERJ5++mnj8f79+6e4n5GRkcYV5bCwMKlYsaLs27dPPRYUFCSBgf/8TIiPj5eYmBjJkyePVKxY0etXJZ966ikVZ8mSRWrXri2O40iPHj3U4/ny5UvUxyxZskj27NlF5J+rsrGxsXL9+nWpXbu20cdly5ZJtmzZjKvggYGB8swzzxj9iI2NlR9++EEeeeQRiYuLk+joaImOjpaYmBhp2bKl7N69W44cOaL6s337dtm9e7dX3wsAQOowsAUAwEtuDFhvDHCTk9wAuGzZsjd9jQMHDkhgYGCibSMiIlLcz1KlSiV6LH/+/HL69GnVTkhIkDfffFPKly8vQUFBEhoaKmFhYbJ161Y5e/Zsil8rNf0JCQmRHDlySGhoaKLH9T6KiMycOVPuvPNOVecaFhYmS5cuNfp44MABKVq0qOTKlcv4W/s927NnjziOI8OHD5ewsDDjv5EjR4qIyMmTJ0VEZMyYMXLmzBmpUKGCVK1aVYYMGSJbt269vTcCAJBq3BUZAAAvCQkJkaJFi950gLN161YpXry45M2b13g8Z86cvuyektydkh3HUfHYsWNl+PDh0r17d3nppZekQIECEhgYKM8++6wkJCT4vD8p6eOsWbMkKipK2rVrJ0OGDJFChQpJlixZZNy4cbJ3795b7seNf9fgwYOlZcuWSW5zYzDcuHFj2bt3ryxatEiWL18u77//vrz55psybdo04wo0ACBtMLAFAMCLWrduLe+9956sXbtW3dlY9+OPP8r+/fuld+/eqXr+0qVLS0JCgvz9999Svnx59fiePXtS3eekzJ07V5o1ayYffPCB8fiZM2cSXUlNL3PnzpVy5crJ/PnzjTtH37i6ekPp0qVl5cqVcvHiReOqrf2e3Vh+KVu2bNKiRYubvn6BAgWkW7du0q1bNzl//rw0btxYRo0axcAWANIBU5EBAPCiIUOGSM6cOaV3794SExNj5GJjY6VPnz6SK1cuGTJkSKqe/8aVxHfeecd4fPLkyanrcDKyZMliXB0VEZkzZ46qMc0IblzV1fv5yy+/yPr1643tWrZsKdeuXZP33ntPPZaQkCBTpkwxtitUqJA0bdpUpk+fLseOHUv0eqdOnVKxvW/z5MkjERERcuXKldT/gwAAqcYVWwAAvKh8+fIyc+ZM6dy5s1StWlV69OghZcuWlf3798sHH3wg0dHRMnv2bAkPD0/V89eqVUsefvhhmTRpksTExEi9evVk9erVsmvXLhGRRGveplbr1q1lzJgx0q1bN6lfv75s27ZNPv30U3VVMyNo3bq1zJ8/X9q3by+tWrWSv//+W6ZNmyaRkZFy/vx5tV27du2kbt26MmjQINmzZ49UqlRJFi9eLLGxsSJivmdTpkyRhg0bStWqVaVnz55Srlw5OXHihKxfv14OHz6s1vGNjIyUpk2bSq1ataRAgQKyYcMGmTt3rvTr1y9t3wQAgIgwsAUAwOs6duwolSpVknHjxqnBbMGCBaVZs2by4osvSpUqVW7r+T/++GMpUqSIzJ49WxYsWCAtWrSQL774QipWrCg5cuTwyr/hxRdflAsXLshnn30mX3zxhdSsWVOWLl0qQ4cO9crze0NUVJQcP35cpk+fLt9++61ERkbKrFmzZM6cObJq1Sq1XZYsWWTp0qUycOBAmTlzpgQGBkr79u1l5MiR0qBBA+M9i4yMlA0bNsjo0aNlxowZEhMTI4UKFZIaNWrIiBEj1HYDBgyQxYsXy/Lly+XKlStSunRpefnll1N9JR4AcHsCHHueEQAAcJ3NmzdLjRo1ZNasWdK5c+f07o4rLFy4UNq3by9r166VBg0apHd3AAC3gRpbAABc5tKlS4kemzRpkgQGBkrjxo3ToUcZn/2excfHy+TJkyVv3rxSs2bNdOoVAMBbmIoMAIDLTJgwQTZu3CjNmjWTrFmzyjfffCPffPON9OrVS0qWLJne3cuQ+vfvL5cuXZK77rpLrly5IvPnz5d169bJ2LFj02yZJQCA7zAVGQAAl1mxYoWMHj1aduzYIefPn5dSpUpJly5dZNiwYZI1K//POimfffaZTJw4Ufbs2SOXL1+WiIgI6du3Lzd7AgA/wcAWAAAAAOBq1NgCAAAAAFyNgS0AAAAAwNVSVIiTkJAgR48eleDgYK8t/I7UcRxH4uLipFixYhIYePv/X4J9m3F4c9+yXzMOjln/xb71X+xb/8V3rX/imPVft7JvUzSwPXr0KHdZzGAOHTokJUqUuO3nYd9mPN7Yt+zXjIdj1n+xb/0X+9Z/8V3rnzhm/VdK9m2K/pdGcHCwVzoE7/HWPmHfZjze2Cfs14yHY9Z/sW/9F/vWf/Fd6584Zv1XSvZJiga2XILPeLy1T9i3GY839gn7NePhmPVf7Fv/xb71X3zX+ieOWf+Vkn3CzaMAAAAAAK7GwBYAAAAA4GoMbAEAAAAArsbAFgAAAADgagxsAQAAAACuxsAWAAAAAOBqDGwBAAAAAK7GwBYAAAAA4GpZ07sDgG7w4MFGO2fOnCq+8847jVyHDh2SfZ6pU6ca7fXr16v4k08+uZ0uAgAAAMhguGILAAAAAHA1BrYAAAAAAFdjKjLS3RdffKFiT9OLbQkJCcnmevfubbRbtGih4tWrVxu5gwcPpvg1kbFUqFBBxTt37jRyAwcOVPHkyZPTrE/4n9y5cxvt1157TcX2Mbpx40aj3bFjRxUfOHDAB70DAAD+hCu2AAAAAABXY2ALAAAAAHA1BrYAAAAAAFejxhZpTq+pFUl5Xa1dQ/ntt9+quFy5ckauTZs2Rjs8PFzFnTt3NnLjxo1L0esj46lRo4aK7Zrrw4cPp3V3YClatKjR7tmzp4rt/VWrVi2j3bp1axVPmTLFB73DzdSsWdNoz58/X8VlypTx+evfe++9RvvPP/9U8aFDh3z++rh1+nfv4sWLjVy/fv1UPG3aNCMXHx/v2475sUKFCqn4yy+/NHLr1q1T8bvvvmvk9u/f79N+2UJCQox248aNVbxs2TIjd+3atTTpE/wPV2wBAAAAAK7GwBYAAAAA4GpMRUaaqF27torbt2+f7Hbbt2832g8++KCKo6Ojjdz58+dVnD17diP3888/G+1q1aqpuGDBginoMdygevXqKr5w4YKRW7BgQRr3BiIiYWFhKp45c2Y69gS3q2XLlkY7KCgoTV/fLinp3r27ih977LE07QuSZn+fvvPOO8lu+/bbb6v4ww8/NHKXLl3ybsf8WP78+Y22/rvJnu574sQJFaf11GMRsz/2km76d4VdirJnzx7fdswP5M2bV8V2SV2VKlVUrC93KeL/07y5YgsAAAAAcDUGtgAAAAAAV2NgCwAAAABwtXStsbWXedGXgjh69KiRu3z5soo//fRTI3f8+HEVMy8/Y9KX/QgICDByen2IXdN17NixFD3/oEGDjHZkZGSy2y5dujRFz4mMR68bETGXj/jkk0/SujsQkQEDBhjtdu3aqbhu3bqpfl59KYjAQPP/wW7ZskXFa9asSfVrILGsWf/3s+CBBx5Ix54krsn797//reLcuXMbObvGHmlDP05FREqUKJHstrNnz1ax/psONxcaGqpie8nEAgUKqNiuce7fv79vO3YT//nPf1RctmxZI9e7d28V89v95uylKl955RUVlyxZMtm/02txRURiYmK827EMhiu2AAAAAABXY2ALAAAAAHC1dJ2KPGHCBKNdpkyZFP2dPn1BRCQuLk7F9nIxaeHw4cMqtv9NGzZsSOvuZEhLlixRcUREhJHT919sbGyqnt9e+iFbtmypeh5kbJUqVTLa+nREe3oW0sabb75ptBMSErzyvA899FCSsYjIgQMHVPzoo48aOXv6Km5Ns2bNVHzXXXcZOfv7zdfsZU30EpNcuXIZOaYipw17yadhw4al+G/1chHHcbzWp8ygZs2aKm7atGmy240ZMyYNepO8O+64w2jrZWL2Enx8Z9+cPrV/0qRJRk5fasvT8TR58mSjrZdwiaT+d3dGxRVbAAAAAICrMbAFAAAAALgaA1sAAAAAgKula42tvryPiMidd96p4j///NPIVa5cWcV6rYGIWW9Qr149I3fo0CEVe7odtu369etG+9SpUyrWl66xHTx40GhTY5uYXh93O4YMGaLiChUqeNz2l19+STKGuzz33HNGW/8scaylna+//lrF9lI8qWUvQXD+/HkVly5d2sjpy0b8+uuvRi5Llixe6U9mYS+hpS/JsnfvXiM3duzYNOnTDW3btk3T18PNVa1a1WjXqlUr2W3t31HffPONT/rkjwoVKmS0H3744WS37dGjh4r136ppRa+r/e6775Ldzq6x1e+vgqQNHjxYxfqyTrfCvg/FfffdZ7T1ZYPsetyrV6+m6jXTE1dsAQAAAACuxsAWAAAAAOBq6ToV+fvvv/fY1i1btizZnL4kQPXq1Y2cvvRDnTp1Uty3y5cvG+1du3ap2J4mrU8PsKduwXtat25ttPXb2mfPnt3InTx50mi/8MILKr548aIPegdfsJcAq127ttHWj0uW+/CdJk2aGO2KFSuq2F7eJ6XL/UybNs1oL1++3GifPXtWxc2bNzdynpYY6du3r4qnTp2aor5kZv/5z3+Mtr6Elj1lTZ8e7iv696n9ufPWUlJIPU9TYm32MY2UmzhxotF+4oknVGwvaTZnzpw06VNyGjVqpOLChQsbuRkzZqh41qxZadUl17LLbrp165bstlu3blXxiRMnjFyLFi2S/buQkBCjrU93/vTTT43c8ePHk+9sBsUVWwAAAACAqzGwBQAAAAC4GgNbAAAAAICrpWuNrbecPn1axStXrkx2O081vDej15XoNb0iItu2bVPxF198kerXgGd2faVdV6uz98Pq1at90if4ll1jZ0uPpQ0yC72++fPPPzdyoaGhKXoOe2mvefPmqXj06NFGzlPtu/08vXr1UnFYWJiRmzBhgopz5Mhh5N5++20VX7t2LdnX83cdOnRQ8QMPPGDk9uzZo+L0WEJLr5+2a2pXrVql4jNnzqRRj6Br3Lixx7y+PIinWnh45jiO0daPhaNHjxq5tFiSJWfOnCp+8cUXjdzTTz+tYrvf3bt3923H/Ix9n6Dg4GAV//jjj0ZO/31kf9d16tRJxfb+Cg8PN9pFihRR8aJFi4zc/fffr+LY2FhPXc8wuGILAAAAAHA1BrYAAAAAAFfzi6nIvlCoUCGj/c4776g4MND8/wH6sjNuuVTvFgsXLlTxvffem+x2H3/8sdG2l7CAO1WtWtVjXp92Cu/KmvV/Xw8pnXosYk77f+yxx4xcdHR0qvpiT0UeN26cit944w0jlytXLhXbn4/FixerODMvzdaxY0cV6++XiPldlxbsJb06d+6s4vj4eCP38ssvqzgzTyVPa/Xr108yToq+7NrmzZt91aVMrVWrVkZbX1bJnqKf2iXP7DKgpk2bqrhevXrJ/t3cuXNT9Xr4R1BQkNHWp3a/+eabyf6dvUTpRx99pGL9fC8iUq5cuWSfxy4JSotp7t7GFVsAAAAAgKsxsAUAAAAAuBoDWwAAAACAq1Fjm4xnnnnGaOtLSujLC4mI/PXXX2nSp8ygaNGiRluv57FrD/R6Pb32SkTk/PnzPugd0oJev9OtWzcj9/vvvxvtFStWpEmfkDx7SRh9eYfU1tTejF4rq9dkiojUqVPHJ6/pZiEhIUbbU41camvyUktfuknErOf+888/jZyn5fzgO7dyTKX158dfvfXWW0a7WbNmKi5WrJiR05dgCggIMHIPPvhgql7ffh57GR/dvn37VGwvLYNboy/TY7Nrq/V70HhiL5Xpyc8//2y03fhbmiu2AAAAAABXY2ALAAAAAHA1piJrGjRooOKhQ4cmu127du2M9h9//OGrLmU68+bNM9oFCxZMdttZs2apODMv3eFvWrRooeICBQoYuWXLlhlt+xb38A17iTPdv/71rzTsyT/0aXJ23zz1ddSoUSru0qWL1/uVUdllHMWLF1fx7Nmz07o7hvDw8GRzfLdmDJ6mMnpreRmYNm7caLTvvPNOFVevXt3I3XfffSoeMmSIkTt16pSKZ86cmeLX/+STT4z2li1bkt123bp1Kua32O2xz8f6VHK7JKBSpUoqtpdGbN++vYrz589v5OxjVs/37NnTyOmfgx07dnjqeobBFVsAAAAAgKsxsAUAAAAAuBoDWwAAAACAq1Fjq3nggQdUnC1bNiP3/fffq3j9+vVp1qfMQK8hqFmzZrLbrVq1ymiPHDnSV11COqpWrZqK7SUG5s6dm9bdybT69Omj4oSEhHTsSWJt2rRRcY0aNYyc3le733qNbWYSFxdntDdv3qxivXZPxKxrj42N9Ul/ChUqpOIOHToku93atWt98vrwrGHDhkb78ccfT3bbs2fPGu3Dhw/7pE+Znb7MpL3sld5+/vnnvfJ65cqVM9r6fQ3084eIyODBg73ymhD57rvvjLZ+fNl1tHrNq6flmOzntJcz/eqrr1Rcvnx5IzdgwAAV678JMjKu2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcLVPX2ObMmdNo62uBXb161cjp9ZzXrl3zbcf8nL027Ysvvqhiu7ZZZ9d1nD9/3qv9QvooUqSI0W7UqJGK//rrLyO3YMGCNOkTzDrW9BAWFqbiyMhII6efMzzR13AUybzn7kuXLhltfa3Jhx9+2MgtXbpUxW+88UaqXq9KlSpG267XK1OmjIo91YZltNruzML+jva0NvSKFSt83R2kgxEjRhht/Ti163jt8yxSz76vwSOPPKJi+x4jISEhyT7P5MmTVWzvr8uXLxvt+fPnq3jo0KFGrmXLliq21xzPqGsWc8UWAAAAAOBqDGwBAAAAAK6WqaciDxkyxGjry0YsW7bMyK1bty5N+pQZDBo0yGjXqVMn2W0XLlyoYpb38U9RUVFGW18K5Jtvvknj3iCjGDZsmIrt5Qk82b9/v4q7du1q5A4ePHjb/fIH+rlUX8ZDRKRVq1Yqnj17dqqePzo62mjb041DQ0NT9DwzZsxI1evj9nhagunMmTNGe/r06T7uDdJCx44djfaTTz5ptPUlw2JiYtKkTzCX6rGPS30ZLvu41KeS21OPbS+99JKKK1eubOT05Tjt6en292tGwRVbAAAAAICrMbAFAAAAALgaA1sAAAAAgKtlqhpbvXZIRGT48OFG+9y5cyoeM2ZMmvQpM/r3v/+d4m379eunYpb38U+lS5dONnf69Ok07AnS09dff220K1asmKrn2bFjh4rXrl17W33yVzt37lSxvpyEiEj16tVVHBERkarnt5elsM2cOVPFnTt3TnY7e5ki+E6JEiVUrNfu2Q4fPmy0N2zY4LM+Ie3cf//9HvNfffWVijdt2uTr7iAJer1tUu3U0s+zX3zxhZHTa2ybNWtm5AoUKKBie5mi9MQVWwAAAACAqzGwBQAAAAC4mt9PRS5YsKCK//vf/xq5LFmyGG19KtzPP//s244hRfSpDteuXUv185w9ezbZ58mWLZuKQ0JCkn2OfPnyGe2UTqmOj4832s8//7yKL168mKLn8GetW7dONrdkyZI07Al0+jIwgYHJ/z9QT1PY3n33XaNdrFixZLe1XyMhIeFmXUxSmzZtUvV3+MfmzZuTjL1p3759KdquSpUqRvuPP/7wRXcgIvXr11exp+NdX4IP/sM+j1+4cMFoT5w4MS27g3Ty5ZdfGm19KvKjjz5q5PRSwYxUvskVWwAAAACAqzGwBQAAAAC4GgNbAAAAAICr+V2NrV03u2zZMhWXLVvWyO3du9do28v/IP1t3brVK88zZ84cFR87dszIFS5cWMV2DYEvHD9+XMWvvPKKz18vI2rYsKGKixQpko49QXKmTp2q4gkTJiS7nb4MhIjn2thbqZtN6bbTpk1L8XMiY9Drt/XYRk1t2tHvR2KLjo5W8VtvvZUW3UEa6NOnj4r130EiIidPnjTaLPGTOdjfu/p3f9u2bY3cyJEjVfz5558buV27dvmgdynDFVsAAAAAgKsxsAUAAAAAuJrfTUUODw832rVq1Up2W3u5FntqMnxDX1ZJJPH0Bl/o2LFjqv7u+vXrKvY0NXLx4sVGe8OGDclu++OPP6aqL/6kffv2KrbLB37//XcVr1mzJs36BNP8+fNVPGTIECMXFhbm89c/deqUiv/8808j16tXLxXbpQXI+BzHSTJG+mnZsmWyuYMHD6pYXzoP7qZPRbaPw6VLlyb7d8HBwUY7f/78KtY/K3A/fcm3ESNGGLnXXntNxWPHjjVyXbp0UfGlS5d807lkcMUWAAAAAOBqDGwBAAAAAK7GwBYAAAAA4Gp+UWNbunRpFS9fvjzZ7ew6MXuZCqSNhx56yGg/99xzKs6WLVuKn+eOO+5Q8a0s0/Phhx8a7f379ye77bx581S8c+fOFL8GTLly5TLaDzzwQLLbzp07V8Xx8fE+6xM8O3DggIofe+wxI9euXTsVDxw40Cevry+FNWXKFJ+8BtJHjhw5ks2ldT1WZmV/19r3J9FdvnxZxdeuXfNZn5Bx2N+9nTt3VvH//d//Gbnt27eruGvXrr7tGNLNxx9/bLR79+6tYvt3/ZgxY1TsrWU7U4ortgAAAAAAV2NgCwAAAABwNb+Yiqwv/VCqVKlkt1u9erXRZpmBjGHChAm3/RyPP/64F3oCX7Gnr50+fVrF9lJJb731Vpr0CSlnL7ukt+3yD/183KZNGyOn7+t3333XyAUEBBjtHTt2pK6zyPC6deum4jNnzhi5l156KY17kznZy9fpS9RVqVLFyO3ZsydN+oSM46mnnjLaPXr0UPEHH3xg5DhmMwd9CT4RkRYtWqjYLul7/vnnVaxPY08LXLEFAAAAALgaA1sAAAAAgKsxsAUAAAAAuJora2wbNmxotPv3759OPQGQEnaNbf369dOpJ/C2ZcuWeWwDtt9++03Fb7zxhpFbuXJlWncnU7KXcxk2bJiK7fuPbNy4MU36hLTVr18/FevLs4gkvq/C1KlTVazfI0NE5OrVqz7oHTK6gwcPqvi7774zcg8++KCKIyMjjZyv75/BFVsAAAAAgKsxsAUAAAAAuJorpyI3atTIaOfJkyfZbffu3avi8+fP+6xPAADg5uxloJD+jh49quLu3bunY0+QVtauXavi5s2bp2NP4HYdOnQw2lu2bFFxRESEkWMqMgAAAAAAHjCwBQAAAAC4GgNbAAAAAICrubLG1hN9XreIyN13363i2NjYtO4OAAAAAPilc+fOGe2yZcumU0+4YgsAAAAAcDkGtgAAAAAAV3PlVORx48Z5bAMAAAAAMg+u2AIAAAAAXI2BLQAAAADA1VI0sHUcx9f9wC3y1j5h32Y83tgn7NeMh2PWf7Fv/Rf71n/xXeufOGb9V0r2SYoGtnFxcbfdGXiXt/YJ+zbj8cY+Yb9mPByz/ot967/Yt/6L71r/xDHrv1KyTwKcFAx/ExIS5OjRoxIcHCwBAQFe6RxSx3EciYuLk2LFiklg4O3PJGffZhze3Lfs14yDY9Z/sW/9F/vWf/Fd6584Zv3XrezbFA1sAQAAAADIqLh5FAAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAABwhf3790tAQIDMmDEjvbsCAMhgGNgCAPzajBkzJCAgQP2XNWtWKV68uERFRcmRI0fSu3te984776T7wC+9+7Bq1SoJCAiQuXPnplsfAABpK2t6dwAAgLQwZswYKVu2rFy+fFl+/vlnmTFjhqxdu1b++OMPyZEjR3p3z2veeecdCQ0NlaioqEzdBwBA5sLAFgCQKdx///1Su3ZtERF56qmnJDQ0VMaPHy+LFy+WRx55JJ17lz4uXLgguXPnTu9uAABw25iKDADIlBo1aiQiInv37jUe37lzp3To0EEKFCggOXLkkNq1a8vixYsT/f2ZM2fk//7v/6RMmTISFBQkJUqUkCeffFKio6PVNidPnpQePXpI4cKFJUeOHFKtWjWZOXOm8Tw36kZff/11effddyU8PFyCgoKkTp068ttvvxnbHj9+XLp16yYlSpSQoKAgKVq0qLRt21b2798vIiJlypSR7du3y+rVq9XU66ZNm4rI/6Zkr169Wp5++mkpVKiQlChRQkREoqKipEyZMon+jaNGjZKAgIBEj8+aNUvq1q0ruXLlkvz580vjxo1l+fLlN+3Djfft2WeflZIlS0pQUJBERETI+PHjJSEhIdH7GxUVJSEhIZIvXz7p2rWrnDlzJlFfUurGv2XXrl3yxBNPSEhIiISFhcnw4cPFcRw5dOiQtG3bVvLmzStFihSRiRMnGn9/9epVGTFihNSqVUtCQkIkd+7c0qhRI1m5cmWi14qJiZEuXbpI3rx5Vd+3bNmSZH1wSj5v165dk9GjR0v58uUlR44cUrBgQWnYsKGsWLEi1e8HAPgbrtgCADKlG4PB/Pnzq8e2b98uDRo0kOLFi8vQoUMld+7c8uWXX0q7du1k3rx50r59exEROX/+vDRq1Ej+/PNP6d69u9SsWVOio6Nl8eLFcvjwYQkNDZVLly5J06ZNZc+ePdKvXz8pW7aszJkzR6KiouTMmTMycOBAoz+fffaZxMXFSe/evSUgIEAmTJggDz30kOzbt0+yZcsmIiIPP/ywbN++Xfr37y9lypSRkydPyooVK+TgwYNSpkwZmTRpkvTv31/y5Mkjw4YNExGRwoULG6/z9NNPS1hYmIwYMUIuXLhwy+/b6NGjZdSoUVK/fn0ZM2aMZM+eXX755Rf54Ycf5N577/XYh4sXL0qTJk3kyJEj0rt3bylVqpSsW7dOXnjhBTl27JhMmjRJREQcx5G2bdvK2rVrpU+fPlK5cmVZsGCBdO3a9Zb7a3v00UelcuXK8uqrr8rSpUvl5ZdflgIFCsj06dOlefPmMn78ePn0009l8ODBUqdOHWncuLGIiJw7d07ef/996dSpk/Ts2VPi4uLkgw8+kJYtW8qvv/4q1atXFxGRhIQEadOmjfz666/St29fqVSpkixatCjJvqf08zZq1CgZN26cPPXUU1K3bl05d+6cbNiwQTZt2iT33HPPbb8nAOAXHAAA/NhHH33kiIjz3XffOadOnXIOHTrkzJ071wkLC3OCgoKcQ4cOqW3vvvtup2rVqs7ly5fVYwkJCU79+vWd8uXLq8dGjBjhiIgzf/78RK+XkJDgOI7jTJo0yRERZ9asWSp39epV56677nLy5MnjnDt3znEcx/n7778dEXEKFizoxMbGqm0XLVrkiIizZMkSx3Ec5/Tp046IOK+99prHf+8dd9zhNGnSJNn3oWHDhs7169eNXNeuXZ3SpUsn+puRI0c6+k+F3bt3O4GBgU779u2d+Pj4JP/dnvrw0ksvOblz53Z27dplPD506FAnS5YszsGDBx3HcZyFCxc6IuJMmDBBbXP9+nWnUaNGjog4H330UXL/fMdxHGflypWOiDhz5sxJ9G/p1auX8ZwlSpRwAgICnFdffVU9fvr0aSdnzpxO165djW2vXLlivM7p06edwoULO927d1ePzZs3zxERZ9KkSeqx+Ph4p3nz5on6ntLPW7Vq1ZxWrVp5/DcDQGbHVGQAQKbQokULCQsLk5IlS0qHDh0kd+7csnjxYjUdNzY2Vn744Qd55JFHJC4uTqKjoyU6OlpiYmKkZcuWsnv3bnUX5Xnz5km1atXUFTXdjam7X3/9tRQpUkQ6deqkctmyZZMBAwbI+fPnZfXq1cbfPfroo8bV4xtTpfft2yciIjlz5pTs2bPLqlWr5PTp06l+H3r27ClZsmRJ1d8uXLhQEhISZMSIERIYaP6ESGrKsm3OnDnSqFEjyZ8/v3p/o6OjpUWLFhIfHy9r1qwRkX/eu6xZs0rfvn3V32bJkkX69++fqn7rnnrqKeM5a9euLY7jSI8ePdTj+fLlk4oVK6r3/sa22bNnF5F/rsrGxsbK9evXpXbt2rJp0ya13bJlyyRbtmzSs2dP9VhgYKA888wzRj9u5fOWL18+2b59u+zevfu2//0A4K+YigwAyBSmTJkiFSpUkLNnz8qHH34oa9askaCgIJXfs2ePOI4jw4cPl+HDhyf5HCdPnpTixYvL3r175eGHH/b4egcOHJDy5csnGgBWrlxZ5XWlSpUy2jcGuTcGsUFBQTJ+/HgZNGiQFC5cWOrVqyetW7eWJ598UooUKZKCd+AfZcuWTfG2tr1790pgYKBERkam6u93794tW7dulbCwsCTzJ0+eFJF/3puiRYtKnjx5jHzFihVT9bo6+30OCQmRHDlySGhoaKLHY2JijMdmzpwpEydOlJ07d8q1a9fU4/p7eqPvuXLlMv42IiLCaN/K523MmDHStm1bqVChglSpUkXuu+8+6dKli9x5550p/4cDgJ9jYAsAyBTq1q2r7orcrl07adiwoTz++OPy119/SZ48edTNiwYPHiwtW7ZM8jnswYk3JXcV1XEcFT/77LPSpk0bWbhwoXz77bcyfPhwGTdunPzwww9So0aNFL1Ozpw5Ez2W3NXW+Pj4FD1nSiUkJMg999wjzz33XJL5ChUqePX1kpLU+5yS937WrFkSFRUl7dq1kyFDhkihQoUkS5YsMm7cuEQ3IEuJW/m8NW7cWPbu3SuLFi2S5cuXy/vvvy9vvvmmTJs2zbgCDQCZGQNbAECmc2NA0qxZM3n77bdl6NChUq5cORH5Z7pwixYtPP59eHi4/PHHHx63KV26tGzdulUSEhKMq7Y7d+5U+dQIDw+XQYMGyaBBg2T37t1SvXp1mThxosyaNUtEUjYl2JY/f/4k7zhsX1UODw+XhIQE2bFjh7pZUlKS60N4eLicP3/+pu9v6dKl5fvvv5fz588bV23/+usvj3/nS3PnzpVy5crJ/PnzjX/fyJEjje1Kly4tK1eulIsXLxpXbffs2WNsdyufNxGRAgUKSLdu3aRbt25y/vx5ady4sYwaNYqBLQD8f9TYAgAypaZNm0rdunVl0qRJcvnyZSlUqJA0bdpUpk+fLseOHUu0/alTp1T88MMPy5YtW2TBggWJtrtxle+BBx6Q48ePyxdffKFy169fl8mTJ0uePHmkSZMmt9TfixcvyuXLl43HwsPDJTg4WK5cuaIey5079y0vixMeHi5nz56VrVu3qseOHTuW6N/Xrl07CQwMlDFjxiRanke/uplcHx555BFZv369fPvtt4lyZ86ckevXr4vIP+/d9evXZerUqSofHx8vkydPvqV/lzfduKqr/zt/+eUXWb9+vbFdy5Yt5dq1a/Lee++pxxISEmTKlCnGdrfyebOnROfJk0ciIiKM/Q4AmR1XbAEAmdaQIUOkY8eOMmPGDOnTp49MmTJFGjZsKFWrVpWePXtKuXLl5MSJE7J+/Xo5fPiwbNmyRf3d3LlzpWPHjtK9e3epVauWxMbGyuLFi2XatGlSrVo16dWrl0yfPl2ioqJk48aNUqZMGZk7d6789NNPMmnSJAkODr6lvu7atUvuvvtueeSRRyQyMlKyZs0qCxYskBMnTshjjz2mtqtVq5ZMnTpVXn75ZYmIiJBChQpJ8+bNPT73Y489Js8//7y0b99eBgwYIBcvXpSpU6dKhQoVjBsjRUREyLBhw+Sll16SRo0ayUMPPSRBQUHy22+/SbFixWTcuHEe+zBkyBBZvHixtG7dWqKioqRWrVpy4cIF2bZtm8ydO1f2798voaGh0qZNG2nQoIEMHTpU9u/fL5GRkTJ//nw5e/bsLb1n3tS6dWuZP3++tG/fXlq1aiV///23TJs2TSIjI+X8+fNqu3bt2kndunVl0KBBsmfPHqlUqZIsXrxYYmNjRcS8mp3Sz1tkZKQ0bdpUatWqJQUKFJANGzbI3LlzpV+/fmn7JgBARpZ+N2QGAMD3bixz89tvvyXKxcfHO+Hh4U54eLhaAmfv3r3Ok08+6RQpUsTJli2bU7x4cad169bO3Llzjb+NiYlx+vXr5xQvXtzJnj27U6JECadr165OdHS02ubEiRNOt27dnNDQUCd79uxO1apVEy1Vc2O5n6SW8RERZ+TIkY7jOE50dLTzzDPPOJUqVXJy587thISEOP/617+cL7/80vib48ePO61atXKCg4MdEVHL7nh6HxzHcZYvX+5UqVLFyZ49u1OxYkVn1qxZiZb7ueHDDz90atSo4QQFBTn58+d3mjRp4qxYseKmfXAcx4mLi3NeeOEFJyIiwsmePbsTGhrq1K9f33n99dedq1evGu9vly5dnLx58zohISFOly5dnN9///22l/s5deqUsW3Xrl2d3LlzJ3qOJk2aOHfccYdqJyQkOGPHjnVKly7tBAUFOTVq1HC++uqrJJdKOnXqlPP44487wcHBTkhIiBMVFeX89NNPjog4n3/+ubFtSj5vL7/8slO3bl0nX758Ts6cOZ1KlSo5r7zyivF+AUBmF+A42pwaAAAAeN3ChQulffv2snbtWmnQoEF6dwcA/A4DWwAAAC+6dOmScffp+Ph4uffee2XDhg1y/PjxJO9MDQC4PdTYAgAAeFH//v3l0qVLctddd8mVK1dk/vz5sm7dOhk7diyDWgDwEa7YAgAAeNFnn30mEydOlD179sjly5clIiJC+vbty82eAMCHGNgCAAAAAFyNdWwBAAAAAK7GwBYAAAAA4GopunlUQkKCHD16VIKDg42FxZH2HMeRuLg4KVasmAQG3v7/l2DfZhze3Lfs14yDY9Z/sW/9F/vWf/Fd6584Zv3XrezbFA1sjx49KiVLlvRK5+Adhw4dkhIlStz287BvMx5v7Fv2a8bDMeu/2Lf+i33rv/iu9U8cs/4rJfs2Rf9LIzg42Csdgvd4a5+wbzMeb+wT9mvGwzHrv9i3/ot967/4rvVPHLP+KyX7JEUDWy7BZzze2ifs24zHG/uE/ZrxcMz6L/at/2Lf+i++a/0Tx6z/Ssk+4eZRAAAAAABXY2ALAAAAAHA1BrYAAAAAAFdjYAsAAAAAcDUGtgAAAAAAV2NgCwAAAABwNQa2AAAAAABXY2ALAAAAAHC1rOndAfgnT4soZ8+e3WhXqlRJxcWKFTNyrVu3VnHevHmN3PHjx1V86dIlI/fXX38Z7RUrVqg4OjrayCUkJCTbV7iHp8+c4zhp2BMkJzAwMMlYRCQ+Pt5os89gH9N8JgAAnnDFFgAAAADgagxsAQAAAACuxlRk+ES2bNmMdtmyZVX83HPPGTlP042zZv3fR9SeuujJlStXjPauXbtU3Lt3byP322+/qZhpyRmLPRUxS5YsRrt06dIqtvfrhg0bVLxkyRIjZ09dh2/Yx/Ozzz6r4ipVqhi59evXG+0ZM2ao+PTp017vG26dfj7OmTOnkdOPKXtauU0/rvXnFBEpUKCAiitUqGDkduzYoWL7M3Gz10Ta0D8X4eHhRi4mJkbFJ06cMHJ896YeZTiZk77f2c//wxVbAAAAAICrMbAFAAAAALgaA1sAAAAAgKtRYwufyJEjh9GuW7euihs3bmzk9NrZq1evGjl9aR67hkCvqbLrvYKDg422XrfVpEkTI7d582YV27W5SF927VCePHmMdt++fVWs12qLmJ+PxYsX+6B3SIp+7Pfo0cPI9evXT8X2MduoUSOjrS8L9tZbbxm5y5cv33Y/cXP20mz6/rTr3b/44gsVnz171sjZ9ZP6cW0f0w0bNlRxp06djNwPP/ygYr0GW0TkwoULifoP37OP4+eff17F9vG/adMmFXfr1s3IxcbG+qB3/kP/naT/nhER6dKli4rt30nfffedivV7jYiIXLt2Ldm/Sy37OzsoKEjF9rEeFxenYvu3HzWjidnn44IFC6rYvm/ImTNn0qJLGRJXbAEAAAAArsbAFgAAAADgahlqKjK3LPcf9rQSfUmHkydPGjl9Oow9fWLp0qUqPnDggJHTt61cubKRe/rpp422Pl2qXLlyHnqOjMQ+J+hTb0RE6tSpo+LcuXMbua1bt6rY/jzCe+zlWu677z4V68v7iIjky5dPxfbyXfYSYU888YSK7Sl0+nmBfetd+jH3wAMPGLmBAweqePbs2UZOP49fv37dyHmainzx4kUjp0+p0z8vIiI1atRQ8ZdffmnkmIqcNjxNNRUxS31CQ0ONnL48260s35cZ2e9zWFiYikeOHGnk2rVrp2K9fEvEXFbJ/g2l/662l8u6leWX9L7aZWh6iUnJkiWN3Ndff51kP+2+ZWb696Jdxqd/DuzSEH2qv/396e/vLWcWAAAAAICrMbAFAAAAALgaA1sAAAAAgKv5vMbW03Id1apVM3Lt27dXsX0L+UOHDqnYXurh8OHDKraXGTh69KiKCxcunGxfRMyaj7x58xq5QoUKqTh//vxG7t1331WxpxqGzMReNmfFihUqtm9VX6ZMGRWvXbvWyOlt+3bmek3X7t27jVzPnj2Ntr78j36LeZHEtSVIX57qdWrXrm20y5Ytq2K7RmfdunUqvpV6Idycvo9q1apl5EaPHq1i/bxpu1kdpn5eePXVV41csWLFVDx16lSPz4Nbo9exv/DCC0ZO/47U65xFRM6dO6fim+0D/XvR/q7Yu3eviu0aTb3O0K6pP3XqlMfXhG+EhIQY7VKlSqnYrqPVfx9RE+2Zfc+BO++8U8XNmzc3cvoyMPZx+dNPP6nYfs+9da7Uvw/sc/6gQYNUbNdjb9++XcX2vVcyK3vMpN8/ZsKECUZOv+eA/X36+uuvq7h79+5Gzt/PlVyxBQAAAAC4GgNbAAAAAICr+WQqsn4p3dOt4F955RUjV758eRXbU1j057Qv1evTKTzdQt5elsKeJqy37Ska+q207Uv+9evXV3HHjh2NnL9f8k+O/d7qU8KnT59u5PRp3/Y0YX0pCPs59f1ZsWJFI1e0aFGjrU93mzt3rpGz9yfSl35861PIRUTq1atntHPlyqXiRYsWGTn9Mwfv0veLvbSWPoXYXoJAn65qTz2zpzTqS72UKFHCyA0ZMkTFO3bsMHKrVq1ScWYtBbkV9ndmp06dVKwvzyIiMnnyZBVv27bNyKV2WqO9j/Rzvv2Z0Kdc2r8tkDbsz0vDhg2NdpEiRVRs71t9eRe7pAwm+3euXrpnT8PfvHmzit977z0jp5fx+ep8qJ/nO3ToYOT08iG9L3abkrB/2FPQGzRooGK99Mpmf17uuOMOFetL54mIfPLJJ0Y7JiZGxf7wnckVWwAAAACAqzGwBQAAAAC4GgNbAAAAAICr+aTGVp+jbc/X1pfmWb9+vZHTl4Gxl3bR62mio6ONnF6TY89P1+s37WV67Bo8vdbSru2LiIhQsV7nI2Le3t5engT/0Ouv7FvO6/vaU420Ta/lGTx4sJGza/v0W95v2rQpxa+BtKd/BvTjTkTkrrvuMtr6Ej8LFiwwcvpyULg99nF59913q/j+++83cvpSbfZSLhs3blTxX3/9ZeT0pVxEzJoyvW5XRKR48eIqnjNnjpHr0aOHihcvXmzk/KF+yNv0WmYRkeHDh6vYPld//PHHKrb3rbfonzW7bkx/zfPnz/vk9eGZ/RurV69eyeb12j0Rkfnz56uYY9Ezu45Wv5eL/ft43rx5KtZ/Y4ukzfusL7/2zDPPGDn997J9PmaJn8TspUb170X9HhU2u35Zv/9I165djZxdq/vhhx+qeMuWLUbOjccpV2wBAAAAAK7GwBYAAAAA4Go+mYqss6crbd++XcX6kg0iIm+++Wayf6dPibJvE6+37alLette7sde5kWf+mFPefzoo49UXLBgQSOnL+lz9uxZgWeepjbY05z0KcX2lLmRI0eqWJ+2KJJ4mtrzzz+vYqaoZmz6VNY6deoYOXu66tatW1W8f/9+n/YrM7Pf9+eee07FegmJ7fjx40ZbX+7DnvJka968uYoff/xxI6cvQ2OXmHzwwQcqrlWrlpE7ePCgit04xcpb9PPq66+/buT077dPP/3UyB07dsy3HbPY39H6ed1XU6GRmP47KjIy0sjpy7mImJ+tJUuWGDm7jAwm/X3Wl7+02/ZvGH36cVqc1+xSr4EDB6rYXmpRLxeaNWuWkeO32D/0/a5PIRYRKVeunIrtJZH03z+zZ882cvpx2axZMyPXvn17ox0eHq7iAQMGGLl9+/ap2C3fmVyxBQAAAAC4GgNbAAAAAICrMbAFAAAAALiaz2tsbfqcent+/a5du1Rs18p6WkIopW5Wk6M/r73cj177aT+PPrfdXh4BN6cvBWTXQevLinTv3j3ZnL0cyapVq4y2XieAjE2v0YmKijJydn29vp/tHG6Pfizec889Rq569eoqtuut9KUoPvnkEyOn19zp9yYQSXwMnz59OslYROTZZ59Vsb70j4hISEiIit944w0j17lzZxVnps+L/X1atWpVFXfq1MnI6cvnTZgwwcjZNV6+6Ju+xJd9zwW9xtYXfUHS9GNcXw5KxLwngojI1atXVTxlyhQj55YavfSin3PtZXP0Glv73HXvvfeqeNu2bUZOX9bSrln39LvaPi7187Nd/9uxY8dkX+O9995TMffBSJr+3tvfg3pNu/3e6kscLl261MitWLFCxfo5XSTxPSsqVKigYvv74LXXXlOxW+5rwBVbAAAAAICrMbAFAAAAALhamk9F9sQb0429pX79+kY7T548KrZvWb98+XIVMz3q9thTz55++mkVN2jQwMjpU6DsqRb60lEi3FbeTfTlY/RppSIicXFxRnv9+vUq5tjzLv34spcHCAoKUrF9rtaX8Zk6daqRi4mJUbFegpCUP/74Q8UHDhwwcnrJx+jRo41caGioiuvVq2fk9Gmu+vP7O7u0ZvLkySrW96WIyKRJk1R88uRJI+eL72X7nF+5cmUV2/3Wz+Pp/RshM9F//zRp0sTI2VNWjxw5ouI9e/b4tmN+Rv++q1GjhpHTpynrS1OKiNx3330qzp49u5HTp+///vvvRu7nn39Wsb4sj0ji30yFChVSsX6OEDHLhw4dOmTk5syZo2J9mjqSZp+P9dIe+3fuDz/8oGJ7XKJ/luxyAfu3Uo4cOVRsT4V24+8qrtgCAAAAAFyNgS0AAAAAwNUY2AIAAAAAXC1D1dimNbs2RL+FeZcuXYycPu/822+/NXK7d+/2Qe8yD30/2LeRr1Spkort2gO9xkq/tbmIyK+//prstshY7OOwRIkSKtbrekTMmiARs/aSfexdnuq99H127tw5I/fqq6+q2K77uZV9pNf22K+h39fArv/VlwErWLCgkWvXrp2Kd+zYYeRuVvPrNvo+su9PULduXRXbSynp9be+qq/S+6bX1IuI9O3bV8V2veCxY8d80h94VqpUKRXnzZvXyNmfEX0JEr0+EDenv18rV640cvq5zD5m9frbpk2bGrl8+fKp+MEHHzRy+hIx9u9au77+ySefVHHz5s2NnP4ZePfdd40cv49vTj8felrqVF+6ScS814S97J6+3+3v4cOHDyfbF3spKbvm1g3c12MAAAAAADQMbAEAAAAArpappyLrt7gWEXnuuedUrF/GFzEv5Y8dO9bIsZTM7cmVK5eKGzdubOT0pSDsKRIbN25U8cCBA40cU6Dcw556U6xYMRXbt7ffvHmz0T5z5oyvupXp2FOOypUrp2JP0w/t6eH6FDpvTQ+3n0c/F9jnBZ2nz5Y+fU/E/5ai0KfxPvroo0ZOP6/u3LnTyF25csW3HRPzu7dt27ZGTi9FuH79upFbs2aNijnHp52OHTuq2D5P2L9/9H3kb9P7fU3/TOu/R0VExo8fr2J72qm+3I5dzqUf6/b3pX4et89/tWvXNtp6yYd97ly1apWK7aWAKBG6NYULFzba+u9je5mn8PBwFRcpUsTI6Uvb2eUe9jJc+ninTp06yfbHXsopo+KKLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVMlWNrV0bUrFiRaOtLwVh1wUsW7ZMxfoSI7h19n7Qaze6du1q5PQlfuxblPfs2TPZHHUd7mF/HurXr69ie7kWe2kJuwYPqWfvh3r16qnYXmpLf98/+eQTI6cvQeAtdk2XpyWh9Lpa+/Nx4sSJJLfzR/r+DA4ONnL6PrLf2xYtWqj4q6++MnJ6DaBdP+lpyQp72R79GNfP43Zf165da+T0ZZ449n1Hr8sUMWts7X1r3wdBv/cFbo1+TNn3Djhy5Eiyf3fw4EEV//LLL8luZ5/jc+bMqeKwsDAj16dPH6Ot5+Pi4oycXg/s6Z4HSJr+e9WuY9XP1fp5U8SsnbXfd71m2q6ftmu09eexa6vHjBmj4n79+iXbt4yEK7YAAAAAAFdjYAsAAAAAcDUGtgAAAAAAV8tUNbZ23UijRo2Mtl5vYNdsDhs2TMVpsc6fPytZsqTR1ufwh4aGGrlTp06peOrUqUbu6NGjyb6GXQdEzW3GZa/Ndvfdd6vY3m87duww2qyT6D123Y2+H+ycflz+/vvvRs5b+0R/TX39WRGRzp07q9het1EXExNjtH/66ScV+/v64/q/b/369UZOr58uVaqUkdPrqOz3XV/z1q7p0veDfUzbtfEtW7ZM8u9EzJrN7777zsgdP3482eeE99h16/ZnRGevianXsSNjsc/Neo2k/fu4UqVKRluvxf/vf/9r5PTvAH5r3R577KGfuyMjI41c1apVVWx/n50+fVrF9jFq/3bW16m3PwcdOnRQ8R9//GHk3nzzTRVnpN9iXLEFAAAAALgaA1sAAAAAgKv5/VRkfUpqjhw5jNxTTz1ltPWpFosWLTJy9vQA3Bp9uYdnnnnGyOnT3exp3voyS0uWLDFy+rb21ONbWcrDG1NnmPqcevbt5fPkyaPi8+fPGzl7Sg3vs/fkypXLaOfPn1/F9udbn+KrT3m6HfayMxUqVFDxv//9byPXpk0bFYeEhBg5fWkDe7mYDRs2qDgjTZ3yBX05nEmTJhk5/TiyS3L09/3RRx81cvpSPPb3qb7/7Omof/75p9HWl5ewp77pS4msWLHCyNlLy8A37O9ofbkv+7iZN2+e0WYZJnfSl1ATESlTpozR1r8D7N9ilAV4j700z+TJk1Vsl3/cf//9Ko6NjTVyetnI1q1bjdzu3buNtv5dP3r0aCNXvXp1FdvnhU8//VTFeplIeuOKLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVMlWNrX5rbBGzlkjEnL8+fvx4I+fv9Vi+ptfBPfTQQ0ZOX9bDrpk8cuSIiu36W33f6jVAST2Pp7ofvT7E3s/6a9g1gHqdWJEiRZJ9/gMHDhjtc+fOJfnamYn+XrZu3drI6fvSrqnVPw/wLvuzbx9TOn1ptMqVKxu5s2fPqtiuF9Jfo2jRokauXbt2RnvgwIEqtuu/7OWHdPpnZtSoUcn2LTOxz3/6PST0+xiIiFSsWFHF9913n5HTz+N6LayIWWu9adMmI2fXxj788MNJvp6IWVdv1+ZSU+87+jnZrrvW33d7X86cOdO3HUOauOeee4y2fq8LEXM5mYxUT+nv9Pd63LhxRu79999XsX2O19v297D9u1O/z4G9rOaECRNUXLhwYSOn1+Pa9bfpWWvPFVsAAAAAgKsxsAUAAAAAuJrfT0XWp9e88sorRs5ermDjxo0qPnbsmG87lsnoy/3YU4r1KQv6diLm9Bh7OoW+dIe9Ly9dumS09SnFdk6f+mZPj4yIiFCxfft7fZmaKlWqGDl9ORR7qpZ++3a9L5lpml2pUqVUrC/dImJ+Hr799lsjZ3924D321CF9ibOyZcsaudDQUBVHRUUZuUKFCqnYnrKmHyddu3Y1cvaUVH26s73ckH4u+OOPP4xc3759VczyUEnT3wd7CYktW7aoeNu2bSl6DhFzH9m5wEDz/6GvXr1axYMGDTJy+pRHjve0U65cORXbx6JeQvDbb78ZuVOnTvm2Y/AZvaSjbt26HrfVSw3OnDnjqy7Bop9L9XOjiMjRo0dT9Bw3W45Sf96FCxcaOf07u1u3bkZOLyn55JNPjNxPP/2U7Ov5GldsAQAAAACuxsAWAAAAAOBqDGwBAAAAAK7m9zW2eh2kXUNg3/J6yJAhyeZwe/SajH379hm5kiVLqtiusb3rrrtUXKdOHSOn1wTaNVx2Ha2et3M6e0kf/Tbodp2AntPrAe3nCQsLM3KelirxV3aNR4cOHVRcvHhxI6cvybJkyRIjl563kPd39nHx2WefqbhmzZpGTl8Konnz5kZOX1bNPp71+lt96RiRxMeFXten18GLmDWaI0aMMHLbt29XMefx23Mry9x5qqOyn0f/rNnnhgsXLqiY49137PddX7rDPjb1/bdgwQIjxzHmXvny5VPxnXfeaeTsY1Y/r9p1+cjYblbjqp9n9VpqEZHp06eruEWLFkZOr8V/4oknjJx+D5y0/rxwxRYAAAAA4GoMbAEAAAAAruZ3U5HtqW9TpkxRsb0kjL0UhH57aniXvmzDsmXLjFytWrVUrC8jImJO97X3rafpFfpUSZun6XX28hL69Ej7Fvf688TFxRm56OhoFf/4449GTl+qJLMsP2JP1X7ooYdUbO9XfcrT1q1bjVxmeb/Sg31cfPXVVyrWb+svItK4cWMVFyhQwMjp09vs6cX2VH+dvZzXoUOHVPz5558bOf28bi83wtTIjE/fR/Z0Y/1zaH9e7M8IUs+ebqwf03Zpjz498bvvvvNtx+Az9vTz8uXLq9j+7WUfa5s2bVKx/fmA/7B/B+zevVvF9lJAAwYMUHGNGjWMXGRkpIr1JeREfP8dzacTAAAAAOBqDGwBAAAAAK7GwBYAAAAA4Gp+UWOr1w3YS8JERESo2K4ZmDRpktGmfsd39DqqmTNnGrnff/9dxZ06dTJyeg1I0aJFjZxef5U7d24jZ8/h15cLses09fpPu+46NjZWxZs3bzZyeh3t/v37jdzff/+t4pMnTxo5u443M8ifP7/RLlasmIrt4+77779Xsb3MC9KOvuzSM888Y+RefPFFFT/22GNGTq9vt2ux9PPAuXPnjJy+PICIyMSJE1W8bt06I+dpyS5kPPY5V7/fhV1Hq3+fUy/tO/pSeiLmOdquxTx27JiKjx8/buS474F76PcsERFp1qyZiu37kly7ds1oBwUFqbhEiRJG7q+//lIxx6x/0Wtu7Rrbtm3bqlj/fIiINGnSRMX2Ep/6d/+tLCmXUlyxBQAAAAC4GgNbAAAAAICr+cVUZP0SeO/evZPN6dNpREQWLVrk244hSfb0Un2aoT3lUGdPj7qVKVD23+r06ZL2c+o5e4oNU7BSTp/SLSKyYsUKFevLw4iIfPPNNyr2xTQV3Dr73PnCCy+o2L6Vf8uWLVVsT086cuSIitesWWPk7GVE9Kn+HGv+5fLlyyq2SxHY12nDfp/1Ehn7+1KfSmhPUYV72Evr6WVYdtmIva0+df3nn382cvqSa/p5W4Tj2e30/WeX6q1cuVLF+rR2EbO0oUiRIkbOXh7T27hiCwAAAABwNQa2AAAAAABXY2ALAAAAAHA1V9bY2vUfNWvWVHHDhg2NnL68xNq1a42cXheAjO92ajU8/a2n29NT4+kd9vIsvXr1UrG9b/Q29TkZk14n//777xu5Dz74IFXPyb7OPPRaLb1OS0Tk4MGDKuYz4Ts//PCD0R47dqyKK1eubOT0+5Ho9dFwF3vf/fLLLyq2l3EqWLCg0daXCrJrJvXlu27nXijI2C5cuGC0ly1bpuKqVasaubCwMBXby3H6+jPBFVsAAAAAgKsxsAUAAAAAuJorpyLbl7X79Omj4tDQUCOnT5lbuHChkdOnKQNIO56mf8PdmHqGmzl79qyKn332WSOnly2wtIzv2Mssvfbaa+nUE6QV+zfvt99+q2J9OR8RkTvuuMNo79q1S8V2GZ/+WaJ8y3/Z+3br1q0q3rZtm5HbvHmzivXykqSex9u4YgsAAAAAcDUGtgAAAAAAV2NgCwAAAABwNdfU2Oq3EC9UqJCRK168uIrPnDlj5Pbu3aviHTt2GDlqAQAASFt6Hfbhw4fTsSdA5qXX3NrL/dhtwL5/hj7emjNnjpHbvXu3iqOjo33aLxtXbAEAAAAArsbAFgAAAADgaq6Ziqw7dOiQ0X7kkUdUnCdPHiN38eJFFZ8+fdq3HQMAAAAAP2JPRdansm/YsMHIBQUFqVgvJU3qebyNK7YAAAAAAFdjYAsAAAAAcLUUTUX29WXjW+2D3R/97sb2nY71tv13GeHflVre6rub3wN/5Y19wn7NeDhm/Rf71n+xb/0X37X+iWM2bXh6fzyN2Xz1mjekaGAbFxd3253xJn1et4hZO5tZ6mjj4uIkJCTEK8+DjMUb+5b9mvFwzPov9q3/Yt/6L75r/RPHbNq4evVqsrkrV6745DVTsm8DnBQMfxMSEuTo0aMSHBycqAgYactxHImLi5NixYpJYODtzyRn32Yc3ty37NeMg2PWf7Fv/Rf71n/xXeufOGb9163s2xQNbAEAAAAAyKi4eRQAAAAAwNUY2AIAAAAAXI2BLQAAAADA1RjYAgAAAABcjYEtAAAAAMDVGNgCAAAAAFyNgS0AAAAAwNX+H5F4mudpopmgAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1200x400 with 16 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "best_model = HVAE(img_size=g_opt.img_size, latent_dim1=g_opt.latent_dim1, latent_dim2=g_opt.latent_dim2)\n",
    "best_model.load_state_dict(torch.load(f\"{g_opt.save_dir}/best.pth\")['model_state_dict'])\n",
    "best_model.to(g_opt.device)\n",
    "best_model.eval()\n",
    "\n",
    "images, labels = next(iter(val_loader))\n",
    "images = images.to(g_opt.device)\n",
    "\n",
    "# 使用模型进行重建（不计算梯度）\n",
    "with torch.no_grad():\n",
    "    output = model(images)\n",
    "    reconstructed_images = output[0]\n",
    "\n",
    "images = images.cpu()\n",
    "reconstructed_images = reconstructed_images.cpu()\n",
    "\n",
    "def show_reconstructions(original, reconstructed, n=8):\n",
    "    plt.figure(figsize=(12, 4))\n",
    "    for i in range(n):\n",
    "        # 显示原始图片\n",
    "        ax = plt.subplot(2, n, i + 1)\n",
    "        plt.imshow(original[i].squeeze(), cmap='gray')\n",
    "        ax.get_xaxis().set_visible(False)\n",
    "        ax.get_yaxis().set_visible(False)\n",
    "        if i == n // 2:\n",
    "            ax.set_title('Original Images')\n",
    "        \n",
    "        # 显示重建图片\n",
    "        ax = plt.subplot(2, n, i + 1 + n)\n",
    "        plt.imshow(reconstructed[i].squeeze(), cmap='gray')\n",
    "        ax.get_xaxis().set_visible(False)\n",
    "        ax.get_yaxis().set_visible(False)\n",
    "        if i == n // 2:\n",
    "            ax.set_title('Reconstructed Images')\n",
    "    plt.show()\n",
    "\n",
    "# 展示前8张图片的重建效果\n",
    "print(\"测试集图像重建效果：\")\n",
    "show_reconstructions(images, reconstructed_images, n=8)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
